% evaluates moment conditions of the form E[h(x,r)] = p 
% aJ is a vector of zeros and ones that selects moments to be evaluated
% set aJ = 0 to use S matrix as covar of pricing errors in case where no
% parameters are estimated (in this case the function always uses all moments) 

function [Jopt,Jfix,S,bvar,bvaropt,pJopt,pJfix,gt,gvar] = GMMcrit(gtfcn,b,x,r,p,W,h,Sfcn,aJ,wt)


t = size(r,1);
if nargin > 7 && size(aJ,1) > 1    % size(aJ,1) == 1 : no SDF parameters estimated
    fast = 0; 
else
    fast = 1; 
end

if nargin < 7
    h = 0; 
end 

[gt, dg] = gtfcn(b,x,r,p,h,fast); 
g = mean(gt);
Jopt = t*g*W*g';    %only a J-stat if inv(S) fed in as W in second stage

if nargin > 6    % use top part only for search of optimum
    
   nq = size(r,2); 
   nf = size(x,2); 

   if nargin < 7
      h = ones(size(r)); 
   end

   if nargin < 9
      aJ = ones(size(g))'; 
   end
   
   if nargin < 10
      wt = 1.5; 
   end 

   if nargin < 8 || isempty(Sfcn)
      S = cov(gt); 
   else 
       
      S = Sfcn( gt-repmat(g,size(gt,1),1) ,h,wt); 
      %S = Sfcn(gt,h); 
   end

   if sum(aJ) == 0 
      covg = (eye(nq) - (dg'/(dg*W*dg'))*dg*W)*S*(eye(nq) - (dg'/(dg*W*dg'))*dg*W)';
      Jfix = t*g*pinv(covg)*g';
   elseif size(aJ,1) == 1
      Jfix = t*(g/S)*g'; 
      gvar = S/t;
   else 
      covg = (eye(nq) - (dg'/(dg*W*dg'))*dg*W)*S*(eye(nq) - (dg'/(dg*W*dg'))*dg*W)'; 
      gind = find(aJ' > 0); 
      subcovg = covg(gind,gind);
      subg = g(gind);
      Jfix = t*subg*pinv(subcovg)*subg';
      gvar = subcovg/t; % variance of pricing error
   end

   bvar = (1/t)*((dg*W*dg')\dg*W*S*W*dg')/(dg*W*dg'); 
   bvaropt = (1/t)*inv(dg*W*dg');

   dgf = min(nq-nf,sum(aJ)); 
 
   pJfix = 1-cdf('chi2',Jfix,dgf); 
   pJopt = 1-cdf('chi2',Jopt,dgf); 

end

